Cjc/vp demo #4320
Cjc/vp demo #4320
Conversation
| -\phi_{x_1x_1} = q_0\int f(x_1,x_2,t)\,\mathrm{d} x_2, | ||
|
|
||
| where :math:`\nabla=(\partial_{x_1},\partial{x_2})`. From now we will | ||
| choose units such that :math:`q_0,m` are absorbed into the definition of |
There was a problem hiding this comment.
Is this a nondimensionalisation? Or simply rescaling f so it has different SI units?
| Each Runge-Kutta stage involves solving for :math:`\phi` before solving | ||
| for :math:`\partial f/\partial t`. Here is the first stage. :: | ||
|
|
||
| # | ||
| fstar.assign(fn) | ||
| phi_solver.solve() | ||
| df_solver.solve() | ||
| f1.assign(fn + df_out) | ||
|
|
||
| The second stage. :: | ||
|
|
||
| # | ||
| fstar.assign(f1) | ||
| phi_solver.solve() | ||
| df_solver.solve() | ||
| f2.assign(3*fn/4 + (f1 + df_out)/4) | ||
|
|
||
| The third stage. :: | ||
|
|
||
| # | ||
| fstar.assign(f2) | ||
| phi_solver.solve() | ||
| df_solver.solve() | ||
| fn.assign(fn/3 + 2*(f2 + df_out)/3) | ||
| t += dt |
There was a problem hiding this comment.
Can this be written as a loop with the coefficient array defined beforehand? I don't think we want to be encouraging people to handcode unrolled Runge-Kutta loops.
There was a problem hiding this comment.
It's a bit fiddly to do that because of the triviality of the first stage, so it will all look a bit cryptic if I do that, and not very pedagogical. This is also how we did it in the original DG advection demo. I think it is harmless.
There was a problem hiding this comment.
Not that I can't use IRKsome because it isn't an IRK.
|
In the future, it could be nice to have a subsequent demo showing how to do this with a single solver, possibly with Irksome too. The whole thing is linear, and if you have a mixed space with (phi, f) then the matrix is lower triangular so you can do exactly the method you have here with a multiplicative fieldsplit to first solve for |
Co-authored-by: Josh Hope-Collins <[email protected]>
Co-authored-by: Josh Hope-Collins <[email protected]>
Co-authored-by: Josh Hope-Collins <[email protected]>
Co-authored-by: Josh Hope-Collins <[email protected]>
Co-authored-by: Josh Hope-Collins <[email protected]>
Co-authored-by: Josh Hope-Collins <[email protected]>
|
Josh, the whole thing isn't linear, because of the a*f appearing in the conservation law. |
Co-authored-by: Josh Hope-Collins <[email protected]>
Co-authored-by: Josh Hope-Collins <[email protected]>
Co-authored-by: Josh Hope-Collins <[email protected]>
Co-authored-by: Josh Hope-Collins <[email protected]>
Co-authored-by: Josh Hope-Collins <[email protected]>
Co-authored-by: Josh Hope-Collins <[email protected]>
Co-authored-by: Josh Hope-Collins <[email protected]>
Co-authored-by: Josh Hope-Collins <[email protected]>
Co-authored-by: Josh Hope-Collins <[email protected]>
Ah yes, of course. I got fooled because each field is linear in itself so you need two linear solves, but the Jacobian of the second solve depends on the solution of the first so just doing forward substitution won't update the lower blocks. |
demos/test_demos_run.py
Outdated
There was a problem hiding this comment.
This file shouldn't have been added here.
demos/vlasov_poisson_1d/vp1d.py.rst
Outdated
| However, the null space also means that the assembled matrix of the | ||
| Poisson problem will be singular, which will prevent us from using a | ||
| direct solver. To deal with this, we will precondition the Poisson problem | ||
| with a version shifted by :math:`\int_{\Omega}\phi\psi\mathrf{d}x`. The |
There was a problem hiding this comment.
| with a version shifted by :math:`\int_{\Omega}\phi\psi\mathrf{d}x`. The | |
| with a version shifted by :math:`\int_{\Omega}\phi\psi\mathrm{d}x`. The |
I think this is your issue. The log says
! Undefined control sequence.
l.17730 ...ifted by \(\int_{\Omega}\phi\psi\mathrf
{d}x\). The
Description
A 1D Vlasov Poisson demo